use strict;
use warnings;
use XML::LibXML;
use Data::Dumper;

my $dir="/home/david/Desktop/busquedas/";

#foreach my $general_index ($xc->findnodes('//mzml:indexedmzML/mzml:indexList/mzml:index')) {
#		my $index_name=$general_index->findvalue('@name');
#		if($index_name eq 'spectrum'){
#			foreach my $index ($general_index->getElementsByTagName('offset')) {
#				my $id_ref=$index->findvalue('@idRef');
#				my $index_location=$index->textContent;
#				$id_ref=~/scan=(\d+)/;
#				my $scan= $1;
#				$indexes{$specs}=$index_location;
#				$specs++;
#			}
#		}
#	}

my %info;
my %prots;



my $file1=$dir."omssa/RGR_2012_331_Fr13_121107_04.xml";
my $file2=$dir."comet/RGR_2012_331_Fr13_121107_04.mgf.pep.xml";
my $file3=$dir."myrimatch/RGR_2012_331_Fr13_121107_04salida.pepXML";
my $file4=$dir."xtandem/RGR_2012_331_Fr13_121107_04.xml";
parse_info($file1,"omssa");
#parse_info($file3,"myri");
#parse_tandem($file4);
#parse_info($file2,"comet");




my $mgf=$dir."spectra/RGR_2012_331_Fr13_121107_04.mgf";
open my $mgf_file,"<",$mgf;
my $specs=0;
while(<$mgf_file>){
	if (/^TITLE\=(.+)\n/){
		$specs++;
		$info{$specs}{name}=$1 if defined $info{$specs};	
	}
}
open my $salida,">","/home/david/Desktop/busquedas/salida.txt";
foreach my $id(keys %info){
	#print $info{$id}{name},"\t";
	my ($xtandem_output,$omssa_output,$comet_output,$myri_output);
	if(defined $info{$id}{xtandem}{seq}){
		$xtandem_output=$info{$id}{xtandem}{seq}."\t".$info{$id}{xtandem}{expect};
	} else{	$xtandem_output="\t";}
	if(defined $info{$id}{omssa}{seq}){
		$omssa_output=$info{$id}{omssa}{seq}."\t".$info{$id}{omssa}{expect};
	} else{	$omssa_output="\t";}
	if(defined $info{$id}{comet}{seq}){
		$comet_output=$info{$id}{comet}{seq}."\t".$info{$id}{comet}{expect};
	} else{	$comet_output="\t";}
	if(defined $info{$id}{myri}{seq}){
		$myri_output=$info{$id}{myri}{seq}."\t".$info{$id}{myri}{mvh};
	} else{	$myri_output="\t";}
	print $salida "$info{$id}{name}\t$xtandem_output\t$omssa_output\t$comet_output\t$myri_output\n";

}	


	################		SIGUE AQUIIIIIIIIIIIIIIIIIIIIIIIIIIIII		###########################################
	################		SIGUE AQUIIIIIIIIIIIIIIIIIIIIIIIIIIIII		###########################################
	################		SIGUE AQUIIIIIIIIIIIIIIIIIIIIIIIIIIIII		###########################################
	#print Dumper %info;
	my %omssa;
	foreach my $id(keys %info){
		$omssa{$id}=$info{$id}{omssa}{expect};
	}
	my @eldest = sort { $omssa{$a} <=> $omssa{$b} } keys %omssa;
	
	my $identifications=0;
	my $num_decoys=0;
	foreach my $id(@eldest){
			$identifications++;
			$num_decoys++ if $info{$id}{omssa}{decoy}==1;
			my $fp=2*$num_decoys;
			#my $tp=$identifications-$fp;
			my $fdr=$fp/$identifications;
			print $identifications,"\t";
			print $id,"\t", sprintf("%.3e",$omssa{$id}),"\t";
			print sprintf("%.5f",$fdr),"\t";
			print $info{$id}{omssa}{decoy},"\t";
			print $info{$id}{omssa}{prots},"\t";
			print "\n";
	}
	#print Dumper @eldest;
	
	################		SIGUE AQUIIIIIIIIIIIIIIIIIIIIIIIIIIIII ARRIBA		###########################################
	
#print Dumper %info;
#print Dumper %prots;
#print "fin\n";

sub parse_tandem{
	my $file=shift;
	#open my $in_file, "<",$file;
	my $parser = XML::LibXML->new();
	my $doc    = $parser->parse_file($file);
	foreach my $spectrum ($doc->findnodes('/bioml/ group')) {
		my $id=$spectrum->findvalue('@id');
		my $decoy_flag=0;
		my @prots;
		
		foreach my $prots ($spectrum->findnodes('protein')) {
			my $prot=$prots->findvalue('@label');
			my ($id,$an)=split /\|/,$prot;
			$decoy_flag=1 if $id eq 'dc';
			push @prots,$an;			
		}
		$info{$id}{xtandem}{prots}=join '|',@prots;
		$info{$id}{xtandem}{decoy}=$decoy_flag;
#		if($info{$id}{xtandem}{decoy}==1){
#			print Dumper $info{$id};
#		}
		foreach my $domain ($spectrum->findnodes('protein/peptide/domain')) {
			my $seq=$domain->findvalue('@seq');
			#print "$id esta repe\n" if defined $info{$id}{xtandem};
			$info{$id}{xtandem}{seq}=$seq;
			$prots{$seq}=join '|',@prots;
			my $expect=$domain->findvalue('@expect');
			$info{$id}{xtandem}{expect}=$expect;
			last;
		}
		
	}
	
	
	
}



sub parse_info{
	my ($file,$engine)=@_;
	open my $in_file, "<",$file;
	my $text_pepxml;
	my $flag=0;
	while(<$in_file>){
		$flag=1 if /\<msms_run_summary/;
		next unless $flag==1;
		$text_pepxml.=$_;
		$flag=0 if /\<\/msms_run_summary/;
	}
	my $parser = XML::LibXML->new();
	#my $doc    = $parser->parse_file($file);
	my $doc = $parser->parse_string($text_pepxml);
	my $num_scans=0;
	foreach my $spectrum ($doc->findnodes('/msms_run_summary/spectrum_query')) {
		my $start_scan=	$spectrum->findvalue('@start_scan');
		my @prots;
		my $decoy_flag=0;
		my $peptide;
		#print "$start_scan\n";
		foreach my $hit($spectrum->findnodes('search_result/search_hit')){
			my $hit_rank=$hit->findvalue('@hit_rank');
			next unless $hit_rank==1;
			if($engine eq 'omssa'){
				push @prots,$hit->findvalue('@protein_descr');
				foreach my $prot($hit->findnodes('alternative_protein')){
					push @prots,$prot->findvalue('@protein_descr');
				}
			}
			else{
				push @prots,$hit->findvalue('@protein');
				foreach my $prot($hit->findnodes('alternative_protein')){
					push @prots,$prot->findvalue('@protein');
				}
			}
			
			
			$peptide=$hit->findvalue('@peptide');
			
#			if($engine eq 'comet' or $engine eq 'myri'){
#				next unless defined $info{$start_scan};
#			}
			if($engine eq 'comet'){
				foreach my $expect($hit->findnodes('search_score')){
					$info{$start_scan}{$engine}{expect}=$expect->findvalue('@value') if $expect->findvalue('@name') eq 'expect';
				}
			}
			if($engine eq 'omssa'){
				foreach my $expect($hit->findnodes('search_score')){
					$info{$start_scan}{$engine}{expect}=$expect->findvalue('@value') if $expect->findvalue('@name') eq 'expect';
				}
			}
			if($engine eq 'myri'){
				foreach my $expect($hit->findnodes('search_score')){
					$info{$start_scan}{$engine}{mvh}=$expect->findvalue('@value') if $expect->findvalue('@name') eq 'mvh';
				}
			}
			$info{$start_scan}{$engine}{seq}=$peptide;
		}
		my @prots_full;
		foreach my $p(@prots){
			my ($id,$an)=split /\|/,$p;
			$decoy_flag=1 if $id eq 'dc';
			push @prots_full,$an;		
		}
		$info{$start_scan}{$engine}{prots}=join '|',@prots_full;
		$info{$start_scan}{$engine}{decoy}=$decoy_flag;
		$prots{$peptide}=join '|',@prots_full;	
		$num_scans++;
	}
	print "Num scans $engine: $num_scans\n";
}


